Surface Nucleation in the Freezing of Gold Nanoparticles. 
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We use molecular simulation to calculate the nucleation free energy barrier for the freezing of 
a 456 atom gold cluster over a range of temperatures. The results show that the embryo of the 
solid cluster grows at the vapor-surface interface for all temperatures studied and that the usual 
classical nucleation model, with the embryo growing in the core of the cluster, is unable to predict 
the shape of the free energy barrier. We use a simple partial wetting model that treats the crystal as 
a lens-shaped nucleus at the liquid-vapor interface and find that the line tension plays an important 
role in the freezing of gold nanoparticles. 

PACS numbers: 61.46.Df, 64.60.Qb, 82.60Nh 



Nanoclusters, consisting of only tens to thousands of 
molecules, exhibit a rich variety of structures and phase 
transitions that are very different from their bulk coun- 
terparts. The most energetically stable structure of a 
cluster varies as a function of the number of atoms or 
molecules. For small clusters, the Mackay or anti-Mackay 
icosahedra are usually the most stable structures, with 
the magic numbers corresponding to completed icosahe- 
dral shells being particularly stable. As the clusters be- 
come larger, decahedra and eventually the face-centred- 
cubic (fee) structures become energetically favorable 
However, molecular dynamics (MD) simulations of cooled 
clusters show that nanoclusters generally freeze to a 
metastable state rather than the energetically most fa- 
vored one, suggesting kinetic factors play an important 
role in determining nanoparticle structure . For exam- 
ple, the fee is the most stable structure for gold particles 
with more than 500 atoms but Bartell et al show that 
clusters with more than N — 1000 still predominantly 
freeze to icosahedra 

Surface phenomena arc expected to play an important 
role in the freezing of nanoparticles. Molecular dynamics 
simulations used to study the freezing of a 561 atom gold 
cluster, cooled at a constant rate from above the melt- 
ing temperature, show the formation of the icoshahedral 
structure is initiated by ordering at the surface rather 
than in the core Q. On the other hand, Lennard- Jones 
clusters initially freeze to a core-ordered icosahedron with 
a disordered surface [H ■ This suggests that the actual lo- 
cation (i.e. surface or core) of nucleation will be driven 
by the wetting behavior of the crystal and liquid inter- 
face. Understanding where nucleation takes place within 
the nanoparticle is important because homogeneous and 
heterogeneous (surface) nucleation rates can differ by or- 
ders of magnitude. 

Furthermore, the location of the forming embryo influ- 
ences the nature of the phenomenological models com- 
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monly used to determine the surface free energy densi- 
ties of crystals. Recent simulations of crystal nucleation 
have shown that classical nucleation theory (CNT) and 
the usual capillarity-based models are not always suffi- 
cient to describe the free energy barrier @, 0|- In the 
case of gold, it is the liquid-solid surface free energy den- 
sity that is of interest and its value is generally obtained 
from models that assume the solid embryo forms in the 
core of the nanoparticle so that it is completely wet by 
the liquid [§]. 

In this letter, we investigate the freezing of gold 
nanoparticles by calculating the free energy barrier to nu- 
cleation using Monte Carlo (MC) simulation techniques. 
Previously, Nam et al [§] calculated the free energy of a 
gold cluster with respect to a global order parameter Q 6 , 
which characterises the total degree of crystallinity of the 
cluster, for a range of temperatures and found that the 
barrier between the liquid phase and the icosahedron was 
considerably lower than the corresponding barrier to the 
more stable fee crystal. We have chosen to characterise 
the nucleation barrier in terms of the size of a solid-like 
embryo n in order to gain insight into the molecular de- 
tails of the nucleation process. This also allows us to 
compare our computed barrier heights with the results 
of CNT, which focuses on the thermodynamic work of 
forming an n-sized embryo. 

Our criterion for identifying an n-sized solid-like em- 
bryo in a cluster is a slightly modified version of the cri- 
terion used in studies of crystal nucleation in bulk sys- 
tems [l^ . Following Frenkel, we begin by defining a 13- 
dimensional complex vector with components, 

n b (i) 

qem ^ = 7T77\ Y ^n{hj) , (1) 

where the sum is over all the neighboring atoms, rib(i), 
within a radius of 3. 5 A of atom i. This distance usu- 
ally contains about 12-13 neighbors if the atom is in the 
core of the cluster. Ye m (d, <p) is the 6 th order spheri- 
cal harmonic and f%j is the unit vector pointing from 
particle i to a neighbor j that specifies the elevation 
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and azimuth angles that this bond makes with respect 
to the coordinate system. Eq. Q] characterises the lo- 
cal order surrounding an individual particle but we also 
expect the local order of neighboring atoms in a solid 
embryo to be highly correlated so we consider two neigh- 
bors to be connected if the correlation function, Cy = 

q 6 W ■ qe(j) = Tf m =-6l6m(i) ■ SemW) where * denotes 
the complex conjugate, is above the threshold value of 
0.65. This threshold value was obtained by comparing 
the distribution functions of obtained for tempera- 
tures above and below the freezing temperature for the 
cluster T — 750K[i(3 and selecting the point at which 
they intersect. The two distributions were not totally 
separated and while the intersect varied slightly, depend- 
ing on the two selected temperatures, this generally oc- 
curred in the range Cy = 0.6 — 0.7. For a much larger 
cluster size (N = 3892) the threshold value can be deter- 
mined more accurately, and is 0.68 [13] • To further distin- 
guish between liquid- and solid-like particles we require 
at least half of the neighbors of a solid particle to be con- 
nected. This last point is where our criterion differs from 
the bulk model and the adjustment is required because 
there are many different environments within a cluster, 
and atoms in the core, surface or at a vertex all have 
very different numbers of neighbors. Finally, two solid 
atoms are considered to be in the same embryo if they 
are connected. Similar order parameters based on Eq. [T] 
have been used previously to identify solid-like particles 
in gold nanoclusters 

We study a N = 465 atom gold cluster, using the 
semi-empirical embedded-atom method (EAM) poten- 
tial [ll|, in the N,V,T ensemble, with a cell volume 
V = 1500A 3 and periodic boundary conditions, to be 
consistent with pervious nucleation studies @, Q. To 
calculate the free energy barrier as a function of embryo 
size n, we implement an umbrella sampling scheme to 
ensure accurate sampling for each embryo size, coupled 
with parallel tempering in temperature 12, 13j] , At each 
temperature, we run eight parallel simulations or win- 
dows, each with a parabolic biasing potential w(n max ) — 
0.0005(n max — no) 2 which biases the system to sample 
states where the largest embryo n max in the cluster is 
around n . We choose n — 0, 10, 20, 30 ... 70 and use 
T = 750,730,710,690,680,670 for tempering. The em- 
bryo criterion is computationally expensive to apply so 
we use trajectories that consist of 10 normal MC moves 
for every particle in the cluster sampling the EAM po- 
tential, followed by a test against w(n max ). If the final 
move is rejected, the system is returned to the state at 
the beginning of the trajectory. During the simulation, 
we attempt to swap clusters with those in neighboring 
windows every 10 trajectories and accept swaps accord- 
ing to the usual replica exchange probabilities [13]. We 
also attempt switches in neighboring temperatures (no 
fixed) every 10 trajectories, but these are offset with the 
n a switches. The tempering switches have acceptance 



ratios of about 0.4 and 0.6 respectively. 

The work of forming an n-sized embryo, Win), is re- 
lated to the probability, P n , of observing the fluctuation 
by a Boltzmann weighted equation. If the system is only 
mildly supercooled, so that the appearance of an embryo 
is rare, then the equilibrium number of clusters, N n , is 
approximately equal to P n so that 
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where k is the Boltzmann constant. We measure the en- 
semble average of N n and use Eq. [5] to obtain W(n)/kT 
from each window. The free energies in each window 
differ by an unknown additive constant, so the full free 
energy curve is constructed fitting the curves to a poly- 
nomial in n [1 21 ] - 

Fig. [U shows the free energy barriers for nucleation at 
four temperatures just below the melting temperature. 
These results were obtained as an average of four inde- 
pendent simulations with each sampling 436,000 trajec- 
tories after the system reached equilibrium i.e. a total of 
1.744 x 10 6 trajectories. An estimate of the error was ob- 
tained by dividing the four runs in half to construct eight 
individual free energy curves and calculating the stan- 
dard deviation of these curves from the average. This 
gives us an error of about 0.6fcT. The free energy curve 
shows the expected maximum and this decreases, along 
with the size of the critical nucleus, with increased su- 
percooling. However, we note that for very small embryo 
sizes, W(n)/kT is actually increasing with decreasing T. 
This feature would not be predicted using a simple CNT 
model. 

We are also interested in understanding where the em- 
bryo forms. In particular, we want to know whether it 
forms at the liquid-vapor interface or within the core of 
the cluster. To identify surface atoms, we use a "cone" 
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FIG. 1: Free energy of formation W(n)/kT of an n-sized em- 
bryo for temperatures T = 730K, 71QK, 69074" and 680A". 
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algorithm [l5[ with an apex angle of 120°. For a cluster 
of 456 atoms, approximately 52% of the atoms are at the 
surface. Fig.[5]shows the average number of atoms in the 
largest embryo found on the cluster surface, n max sur f ace 
versus the embryo size n max . The cross symbols repre- 
sent the values for all the temperatures studied and we 
have made these all the same color and symbol for clar- 
ity. The results for each of the individual temperatures 
can be fitted to a linear curve but all the lines are es- 
sentially the same within the scatter of the data so we 
have averaged over all temperatures to obtain the dark 
symbols. Approximately 46% of the atoms in the larger 
embryos and 63% of the atoms in the smaller embryos lie 
at the cluster-vapor interface, clearly suggesting that the 
growing solid phase is partially wet by the liquid. Fig. [3] 
shows a typical large embryo. 




FIG. 2: The average number of surface atoms Umax, surf found 
in the largest embryo n m ax- The dark circles represent the 
values averaged over all temperatures. The lines are linear 
best fits for the larger (red) and smaller (blue) embryos. The 
cross symbols represent the averages for each temperature 
studied and give an indication of the scatter of the data. The 
insert is a magnification of the small embryo region. In color 
online. 




FIG. 3: Left: An isolated 71 atom embryo. Right: The 
same solid embryo (dark atoms) embedded in the liquid (light 
atoms) cluster. In color online. 

The condition of partial wetting places thermodynamic 



limits on the surface free energy densities (ov,) of the 
three phases such that a 13 — 023 < o\i, where the 
subscripts 1, 2 and 3 denote the solid, liquid and va- 
por phases respectively [l6j |. For the EAM potential, 



(T13 = 0.90J/m 2 and 023 = 0.74J/m 2 [Tjj which requires 



012 > 0.16J/m 2 . Bartell et al [8] found that a num- 
ber of thermodynamic theories and empirical relations 
gave estimates of o\2 in the range of 0.11 — 0.16J/m 2 . 
In the same work, the authors used a CNT model that 
assumes complete wetting of the solid embryo by the liq- 
uid (core nucleation) to predict the solid-liquid surface 
tension based on fitting the rate of nucleation obtained 
from molecular dynamics simulations. For a cluster of 
N = 459, at T = 720K, they found a 12 = 0.084J/m 2 
which is well below the wetting threshold. 

To obtain an estimate of CT12 under the conditions of 
partial wetting, we assume the solid embryo grows at a 
planar liquid-vapor interface in the shape of a lens (see 
Fig. 2]). The lens model, without a line tension contri- 
bution, has previously been used to study droplet forma- 
tion at the liquid-vapor interface [l8| and its application 
to crystal nucleation in a small nanoparticle implies a 
significant simplification of the microscopic process. We 
continue to ascribe bulk-like properties, such as chemi- 
cal potential and surface tensions, to a nucleus containing 
just tens of atoms, in line with CNT, and we ignore many 
effects including surface roughness and the appearance of 
crystal facets. Nevertheless, such a model is numerically 
tractable and provides some useful insight into the nucle- 
ation process. Fits to our simulation data using the core 
nucleation CNT model, or the lens model without the 
line tension, were unable to account for the shape of the 
free energy curve. Auer and Frenkel found that the 
line tension, t, and its curvature correction, To, played 
an important role in the heterogeneous freezing of hard 
sphere colloids at a wall. Including r in the mechanical 
equilibrium at the three phase contact line for the lens 
model makes the contact angles size-dependent [l^l so 
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The work needed to form an n-sized embryo can then 
be expressed 



W(n) = nAfi+R 2 



St 

i=2,3 
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where Afj, is the difference in chemical potential between 
the liquid and solid phase. R is obtained numerically un- 




FIG. 4: Schematic cross section of a solid (phase 1) lens nu- 
cleus forming at the liquid (phase 2) - vapor (phase 3) inter- 
face. A12 and A13 are the solid-liquid and solid- vapor interfa- 
cial surface areas respectively and R is the radius of the lens. 
The four arrows originating from the 3-phase contact are the 
force vectors of the surface tensions Oij and the line tension 
t/R. 




°0 10 20 30 40 50 60 70 80 
n 

FIG. 5: Phenomenological model data fits to the calculated 
free energy barrier at T — 710K. The lens model with line 
tension included (solid line) and CNT, assuming core nucle- 
ation (dotted line). 



der the constraint of fixed volume for an n-sized embryo, 
77 d3 a a sin0ii(2 + cos0ii) 

i=2,3 y ' 

where v = 1.7277 x 10 _29 m 3 is the volume per molecule 
in the solid phase. 

Fig. [5] shows a fit of the model to the data at T = 710 
where A/x/fcT = -0.22, cr 2 3 = 0.18J/m 2 , r = -1.17 x 
10 _11 J/to and To = 3.92 x 10~ 21 J were adjustable pa- 
rameters. Due to the nature of the numerical fit, we 
cannot guarantee that we are at the global minimum for 
the data fit and there is still the possibility of finding an 
improved fit with a positive line tension. At small embryo 
sizes, we can expect a greater fraction of particles to be 
in the surface of the lens based on surface-to-volume ra- 
tio arguments, and hence, a greater fraction of solid-like 
particles in the small embryos will appear on the surface 
of the cluster, as compared to the larger embryos. A neg- 
ative r would enhance this effect by stretching the lens 
of small embryos, where line tension is most important. 
We also fit the CNT core nucleation model to our data 
at T — 710K, using /i and a 12 as adjustable parameters 
and assuming a spherical geometry for the embryo. The 
resulting a 12 — 0.085 J/m 2 is the same as that obtained 
from direct measurements of the rate [1]. We see in Fig. [5] 
that the model clearly fails to predict the correct shape of 
the barrier but does obtain a close estimate of the barrier 
height. 
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